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ABSTRACT 

The problem of the determination of the subcritical 
aeroelastic response and flutter instability of nonlinear 
two-dimensional lifting surfaces in an incompressible 
flow-field via Volterra series approach is addressed. 
The related aeroelastic governing equations are based 
upon the inclusion of structural nonlinearities, of the 
linear unsteady aerodynamics and consideration of an 
arbitrary time-dependent external pressure pulse. 
Unsteady aeroelastic nonlinear kernels are 

determined, and based on these, frequency and time 
histories of the subcritical aeroelastic response are 
obtained, and in this context the influence of 
geometric nonlinearities is emphasized. Conclusions 
and results displaying the implications of the 
considered effects are supplied. 

NOMENCLATURE 

a Dimensionless elastic axis position measured 
from the midchord, positive aft 
c Chord length of 2-D lifting surface, 2b 

Chi * Ctf * Kfc , Damping and stiffness coefficients in 

plunging and pitching (i=l,2,3 - linear, 
quadratic, cubic), respectively 
C La Lift-curve slope 

C(k\ F(k\G(k) Theodorsen’s function and its real and 
imaginary counterparts, respectively 
h,£ Plunging displacement and its dimensionless 
counterpart, (h/b), respectively 
h n , H n n -th order Volterra kernel in time, and its 
Laplace transformed counterpart, respectively 
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/ , r a Mass moment of inertia per unit wingspan and 
the dimensionless radius of gyration, 

(l a /mb 1 2 )'' 2 , respectively 

l ,m a Dimensionless aerodynamic lift and moment, 
(i LJb/mUi ) and (m „ b 2 /l a Ul), respectively 
L , M (l Total lift and moment per unit span 
L h , l h Overpressure signature of the N-wave shock 

pulse and its dimensionless counterpart, 

( L b b/mUl ). respectively 

m,fi Airfoil mass per unit length and reduced mass 
ratio, (m/xpb 2 ), respectively 
N Load factor, 1 + h*/g 

P ,p Peak reflected pressure amplitude and its 
dimensionless counterpart, ( P m b/mUl)i 
respectively 

r Shock pulse length factor 

s j , X Laplace transform variable and Laplace 
operator, respectively, Sj =ik j ;i 2 =-l 
S , y Static unbalance about the elastic axis and its 

Ao 

dimensionless counterpart, S a j mb » 
respectively 

t , T q , TTime variables and dimensionless counterpart, 
{Uj/b), respectively 

tp,?p Positive phase duration, measured from the 

time of the arrival of the pulse, and its 

dimensionless value, respectively 
TF Transfer function 

, V Freestream speed and its dimensionless 
counterpart, (U„/bo) a ) 

x(t) Time-dependent external pulse (traveling gusts 
and wake blast waves) 

y(t) Response of the considered degree of freedom 

(pitch a and/or plunge h) 
a Twist angle about the pitch axis 

£ h >£ a Structural damping ratios in plunging 

(c h /2m<o„), and pitching (c a /2I a (o a ), 
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respectively 
p Air density 

Wagner and Kussner’s functions in the time 
domain, respectively 

CO , k Circular and reduced frequencies, {(obfU J), 
respectively 

CO h , co a Uncoupled frequencies in plunging and 
pitching, ( Kjm )" 2 and ( 
respectively 

i 0 Plunging-pitching frequency ratio, (coJco a ) 

Superscript 

( . ) Variables in Laplace transformed space 

( . ), (. ) Derivatives with respect to the time t , and the 
dimensionless timer , respectively 

INTRODUCTION 

I t is a well - known fact that within the linearized 
approach of the aeroelasticity discipline, it is possible 
to obtain the divergence and the flutter instability 
boundaries, and also to get the linearized subcritical 
aeroelastic response of flying vehicle exposed to time- 
dependent external pulses. On the other hand, in 
addition to the above mentioned items, the nonlinear 
approach of the problem can provide important 
information about the nature of the instability boundary, 
i.e. benign or catastrophic one, and on the influence of 
the considered nonlinearities on the subcritical 
aeroelastic response. In other words, such an approach 
gives the possibility of determining in what conditions 
the flutter speed can be exceeded without the occurrence 
of a catastrophic failure of the lifting surfaces, in which 
case the flutter is benign , as well as the conditions in 
which undamped oscillations may appear at velocities 
below the flutter velocity, in which case the flutter is 
catastrophic. Due to the strong implications of various 
nonlinearities on the highly flexible lifting surfaces, 
their related aeroelastic phenomena cannot longer be 
analyzed solely within the standard linearized 
aeroelasticity theory. Aircraft wing structures often 
exhibit nonlinearities, which affect their aeroelastic 
behavior and performance characteristics and flutter 
boundaries. In order to investigate the aeroelastic 
behavior of the aircraft in the subcritical flight speed 
range, and in the vicinity of the flutter boundary, the 
aeroelastic governing equations have to be considered in 
nonlinear form. 

This investigation concerns the time and frequency 
formulations of nonlinear two-dimensional lifting 
surfaces exposed to an incompressible flow field and 
subjected to an external pressure pulse 1 ' 3 . 


Based on Volterra’s functional series approach 4 ’ 8 
important information about the effects of nonlinearities 
on either the aeroelastic response in the subcritical flight 
speed regime, and their implication on the benign or 
catastrophic character of the flutter boundary are 
supplied. 

The advantage of the technique based on Volterra’s 
series and indicial function (Lomax 9 , Bisplinghoff 10 
Marzocca et al. 11 ’ 12 ) consists, among others, on the 
possibility to investigate, within a rigorous theoretical 
basis, the aeroelastic systems featuring a wide class of 
structural nonlinearities. 

First of all, based upon the first order Volterra kernel 
the study of the aeroelastic stability of the systems can 
be carried out. Moreover, this methodology can 
encompass the case of an arbitrary number of degrees of 
freedom and at the same time is conceptually clearer, 
computationally simpler and can provide more accurate 
and realistic results as compared to the conventional 
techniques used in nonlinear aeroelastic systems based 
on perturbation and multiple scale methods. 

Toward the end of determining the nonlinear unsteady 
aeroelastic kernels, the harmonic probing algorithm, 
referred to as the method of growing exponentials 
advanced by Bedrosian and Rice 13 , and the 
multidimensional Laplace transform will be used. 

In addition to the aeroelastic response and determination 
of the flutter instability boundary, Volterra Series will 
be used to study the conditions rendering the flutter 
boundary a benign or a catastrophic one (Librescu 14 ' 15 ). 
Moreover, when the closed-loop dynamic response of 
actively controlled lifting surface is analyzed, also the 
feedback control forces and moments will be included 
(Librescu 15 , Librescu and Gern 16 , Librescu and Na 17 , 
Van Trees 18 , Chua and Ng 19 ). 

The Volterra’ s series approach provides a firm basis of 
nonlinear subcritical aeroelastic response, in the sense 
that it supplies an explicit relationship between the input 
(any kind of time-dependent external pulses, i.e. blast 
load, sonic-boom, gust loads) and its response. 

With the so-called Volterra Kernel identification 
scheme the modeling of an aeroelastic system using this 
approach becomes feasible. However, this methodology 
requires determination for each specific flight 
conditions of the appropriate nonlinear kernel of the 
Volterra’s series. For this reason, in order to define the 
appropriate aerodynamic loads, the recent interest in the 
modeling of unsteady nonlinear aerodynamics by this 
approach has been focused on the identification of 
Volterra’s kernels in the time domain (Silva 20 ' 23 ), and in 
the frequency domain (Marzocca at al. 12 , Tromp and 
Jenkins 24 ). 

A number of fundamental contributions related with 
Volterra’s series, developed by outstanding 
mathematicians (Volterra 4 , Wiener 5 ) and used mainly in 
electrical engineering 6 ’ 8 , are already available. 
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The original studies on functional series by Volterra 4 , 
have been continued in the works by Volterra himself, 
of those of famous physicists and mathematicians as 
Rugh 6 , Schetzen 7 , Boyd 8 . These concepts have been 
used in nonlinear system theory, in general, and in the 
modeling of nonlinear aeroelastic systems, in particular 
(Silva 20 ). Very few applications of this method have 
been done in the aeroelasticity discipline. 

Originally, the method of Volterra series and Volterra 
kernel identification were developed to identify the 
nonlinear behavior in electrical circuits. In the aerospace 
field, the fundamental contributions were brought by 
Silva, who has shown that the method is also applicable 
to aeroelastic systems (aerodynamic reactions and 
forced structural model). 

Silva’s pioneering work 20 23 , in this area has opened a 
very promising way of modeling and approaching 
nonlinear aeroelastic systems. 

BASIC CONCEPTS AND LIMITATIONS 
OF THIS APPROACH 

Having in view the fact that for nonlinear systems the 
superposition principle is not applicable, and having in 
view the different types of responses induced by 
unsteady aerodynamic loads and the external excitation, 
a combination of transfer functions is used. These 
transfer functions for the nonlinear aeroelastic systems 
and the time-histories response in time and frequency 
domains are determined by taking the multi-dimensional 
Laplace transform of the Volterra kernels of the related 
aeroelastic system via a Mathematical code developed 
by these authors 25 . 

Our approach intended to address the subcritical 
response of the nonlinear aeroelastic governing 
equations, is based on its exact representation as an 
infinite sum of multidimensional convolution integrals, 
the first one, (i.e. the linear kernel) being the analogous 
to the linear indicial aeroelastic function. The full 
nonlinear aeroelastic response will be composed of 
additional higher-order contributions. In the frequency 
domain, if the nonlinear function governing a system is 
’smooth’, then for small inputs the system must be 
asymptotically linear 6 . One of the key issues is to 
determine, corresponding to the considered type of 
structural and aerodynamic nonlinearities, the pertinent 
Volterra's kernels. When also the active control is 
implemented the corresponding Volterra’s kernel should 
also be derived. 

THE THEORY 

In an attempt to make the paper as self-contained as 
possible, several elements associated with Volterra’s 
series as applied to aeroelastic system, as well as with 
the indicial functions will be supplied here. 


Indicial Theory and Aerodynamic Loads 

Using the aerodynamic indicial functions 
corresponding to transient aerodynamic reaction to a 
step pulse, the aerodynamic forces and moments 
induced in any maneuver and any flight regime can be 
determined. Aerodynamic forces and moments acting 
on a rapidly maneuvering aircraft are, in general, 
nonlinear functions of the motion variables, their time 
rate of change, and the history of the maneuvering 
(Tobak & Chapman 26 ). However, in this study, the 
linear aerodynamic theory is adopted. 

Once the response of the system to a step change in one 
of the disturbing variables (i.e. the indicial response) is 
known, the indicial method permits the determination of 
the response of a system to an arbitrary schedule of 
disturbances. There is a critical value of the flight speed 
above which the steady motion becomes unstable. In a 
nonlinear aeroelastic system the flutter phenomenon 
corresponds to the instability known as the Hopf 
bifurcation, resulting in finite amplitude oscillations, in 
the case of supercritical Hopf bifurcation , and in 
oscillations with increasing amplitudes, even if the 
system operates before reaching the flutter speed, in the 
case of the subcritical Hopf bifurcation 21 ' 29 . 

We need to mention that a nonlinear indicial theory 30 , 
asserts that the response of a nonlinear system to an 
arbitrary input can be constructed by integrating a 
nonlinear functional, that involves the knowledge of the 
time-dependent input and the kernel response. Whereas, 
within the linear indicial theory the linear kernel or 
linear impulse response can be convolved with the input 
to predict the output of a linear system, the nonlinear 
indicial theory constitutes a generalization of this 
concept. It can also be shown that the traditional 
Volterra-Wiener theory of nonlinear systems constitutes 
a subset of nonlinear indicial theory. It should also be 
mentioned that the nonlinear unsteady aerodynamics 
valid throughout the subsonic incompressible/ 
compressible, transonic and supersonic flight speed 
regimes can be used and determined via the use of 
nonlinear indicial functions 31 in conjunction with the 
Volterra’s series approach. 

Volterra Functional Series Theory 

As it was shown (Rugh 6 , Schetzen^) within Volterra’s 
series approach the full response in the time domain, 
y(t)> of the nonlinear systems with memory can be cast 
as: 

y(0=iy*(/)- <•> 

*=0 

where, y (t) is expressed as: 

>’»(')= fJJJT (2) 

k limes »=1 

By a change of variables, it is possible to express Eq. 
(2) in contracted form as: 


3 

American Institute of Aeronautics and Astronautics 



?*(')= IJJr K (f, ,T 2 -r* )fl x{t - r, )dr, ■ (3) 

k timet ,= * 

It is assumed that x(t) = 0 for r<0 implying that the 
system is causal. 

With this restriction, all the integrals in the subsequent 
discussions are different from zero over the time range 
(0,°o). Restricting the development of Eq. (3) to the 
first three terms one obtains: 

_v(f) = Jy^(r l )x(/-T| ]dz, 

+ JJ >h (r, ■ h )x{t - r, )x(t - r 2 \h x dx i 

+ [f| h, (t, , r 2 , t 3 )x(t - r, )x(t - r 2 )x(t - t, ]dT,dT 2 dt 3 

+ -(4) 

On the other hand, the response of the system can be 
expressed also in the frequency domain. 

The Volterra series is essentially a polynomial 
approximation of the system, extension of Taylor series 
to systems with memory, while Volterra 1 s kernels 
h i (s i ) are a direct extension of the impulse response 
concept of the linear system theory to multiple 
dimensions (Volterra 4 5 , Rugh 6 , Schetzen 7 , Boyd 8 ). 
Consequently, a multidimensional analogue of the 
impulse response can be used to characterize a nonlinear 
system (Silva 20 ' 23 ). 

Having in view that the aeroelastic systems memory is 
not infinite and, at the same time, the time-dependent 
external excitations, such as impulse, gust, blast and 
sonic-boom pressure signatures are non persistent but 
their effect will diminish as time unfolds, it is possible 
to characterize a nonlinear aeroelastic system via 
Volterra series. This fact is reflected in the 
interpretation of the Volterra kernels as higher order 
impulse response functions, i.e. /i(tj „)— > 0 as 


->°o. 

We will use the definition of the nonlinear transfer 
function or higher-order impulse response functions 
namely: 

H k (sj , s 2 , s,, ) = JJ-*- J h k (r P T 2 >* ) 

e -sih e -s 2 r i ... e ~ x ' ,iTHi dr l dT 2 ••• dx k ,( 5 ) 
as well as of its inverted counterpart: 


2ju Ja .- i 00 Jo 2 - f °° 


e^'e* 1 * 2 -••e s ’ ,rm ds l ds 2 --‘ds n - (6) 


Once the Volterra’s kernels are known the response of 
the nonlinear aeroelastic system can fully be identified. 
As demonstrate in the Schetzen works 7 , without loss of 
generality, the kernels will be taken as symmetric. 

If we focalize the attention on the linear system, the 


Laplace transform X of the first term of Eq. (4) yields 
the familiar Laplace domain expression y( s )~ H(s)x(s) 


where y(s\h(s\x(s) are the Laplace transforms of 
x(t), respectively, and H(s) is the transfer 
function of the system; either the first transfer function 
or the first kernel in time h(z) encode all the 
information about the aeroelastic system, that is, of 
course, exact only for the linear system. Moreover, as is 
well known, if the system is linear, i.e. superposition 
principle holds valid, and is time invariant, the external 
load is uniquely related to the response by a convolution 
integral. With the use of functional series, i.e. the 
Volterra series, this functional representation can be 
extended to nonlinear systems. The comparison between 
the prediction of the linear aeroelastic responses of 2-D 
lifting surface in incompressible flow field based on the 
Volterra’s series approach (using Theodorsen’s 
function) and on the exact solution, based on 
convolution integrals (using Wagner’s function) is 
presented in Fig. 1 . 



» 2 4 6 $ 10 


Time [sec] 

Fig. 1 Aeroelastic response time-history of lifting 
surface in incompressible flight speed regime to Dirac 
delta impulse, as represented in inset. Comparison of 
prediction of responses based on the first Volterra 
kernel and the exact solution. 

The excellent agreement of these two predictions shown 
here, assess both the accuracy of the aeroelastic model 
and also the power of the methodology that combines 
Volterra’s series and indicial function. 

MATHEMATICAL FORMULATION 


General theory for 2-D lifting surfaces 
including structural nonlinearities 

The aeroelastic governing equation of motion for 1 
and 2 DOF including structural nonlinearities that 
include the damping and the stiffnesses can be analyzed 
in the following way. Two systems will be analyzed 
here: a 1 DOF lifting surface (i.e. plunging only) and a 2 
DOF lifting surface featuring structural and 
aerodynamic coupling in plunging h and pitching a. As 
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previously mentioned, the unsteady aerodynamic is 
considered linear. A harmonic time dependent external 
concentrated load is also applied. This configuration for 
example, can be considered to correspond to an engine 
mounted on an aircraft wing. As a result, a harmonic 
type loading due to the engine oscillations has to impact 
the motion of the wing. 

As a characteristic of this approach, the transfer 
functions of the system would exist and be the same for 
any excitations 32 ' \ (namely for random, sine, impulse). 
This is due to the fact that transfer functions are a 
characteristic of the system itself and are independent of 
the input to the system. 

As a reminder, the validity of this method is based on 
the use of continuous polynomial type nonlinearities. 

For nonlinear ordinary differential systems, there are in 
general, an infinite number of Volterra kernels. In 
practice, one can handle only a finite number of terms in 
the series, which leads to the problem of truncation 
accuracy. However, Wiener suggest that the first terms 
of the series may be sufficient to represent the output of 
a nonlinear system if the nonlinearities are not too 
strong. 

The use of the multidimensional Laplace transform as a 
function of several variables is a tool useful in 
stationary nonlinear system theory. The multivariable 
convolutions can be represented in terms of products of 
Laplace transforms. 

It is well known that the nonlinear aeroelastic systems 
cannot be described by a simple transfer function for 
two main reasons: a) the response has different trends as 
compared to the unsteady aerodynamic loads and the 
external excitation and, b) in the nonlinear case the 
superposition principle is not applicable. It is also well 
known that any time-dependent external excitation, i.e. 
periodic or otherwise, can be represented, to an arbitrary 
degree of accuracy, by a sum of sinusoidal waves 31 . In 
this context, if the external load is expressed in term of 
multiple sinusoidal form (for example traveling gust 
loads) this is easily convertible in the exponential form, 
i.e.: 

u(t) = Acos(a) A t) + Bcos(co B t) 

<=> u(t) = ^-{e SAt + e ~ Sfi ) . ( 7 ) 

For clarity of exposition, it is convenient to adopt this 
approach for a system with one degree of freedom 
(l DOF), These results have more general bearing and 
can be extended for systems with multi-degree of 
freedom (MDOF). In fact, by using the classical 
approach of the one dimensional frequency response 
function, it is possible to derive an analytical form of 
the multi - dimensional frequency response 
characteristics of nonlinear systems. The systems based 
on 1 DOF (plunging h) and 2 DOF, pitching a and 
plunging /i, will be considered in the next sections. 


Pure Plunging Airfoil 

The nonlinear equation of motion of an airfoil 
featuring plunging motion can be expressed as: 

7=1 

where k hh c hh m are the stiffness, damping and mass 
parameters, respectively. 

In Eq. (8) the related unsteady aerodynamic lift is 
represented as a function of the plunging degree of 
freedom h, only: 

Kb) = -C La pUl\[j{x-z»ydx u -\pC La Ulh’ ■ (9) 

The non-circulatory component present in Eq. (9) has 
been represented in term of convolution integral of the 
indicial Wagner’s function. 

In order to explain how this methodology works, let us 
to determine, in terms of Volterra series, how a system 
responds to a harmonic or periodic time-dependent load. 
Let consider a periodic external excitation in the form: 

<l0) 

As is well known, the information acquired by the case 
of the response to a harmonically time-dependent load 
can be used to obtain the response to any time- 
dependent excitation. In fact, considering the case of a 
concentrated load arbitrarily located in the x, y plane of 
the wing, we have: 

u{ x, yj)=A8{x- x„, y - y 0 >'• , (11) 

where <$(.), x 0 , y 0 , A, co denote Dirac’s distribution, 

location of the load, its amplitude and excitation 
frequency, respectively. Once determined the transfer 
function (labeled as TF) corresponding to a given 
excitation frequency, its counterpart in the time domain 

can be obtained as the inverse Laplace transform £~ ] : 

TF(x, \\t) = £ -'{TF(x, y,.v)}= — y,s)e"ds 

2 

In addition to the direct role in the determining of the 
response, the transfer function TF has then the role in 
determining the response to arbitrary time-dependent 
external excitations. 

The general procedure to identify the aeroelastic kernels 
of various order is to consider a general input in 
the form of Eq. (10) and to equate coefficients of 
X x X 2 -'-X tl e { ' ,+S2+ "+***. As an example, the first 

aeroelastic Volterra’ s kernels that describes the linear 
system in the aeroelastic governing equations, obtained 
by neglecting the nonlinear terms, is obtained by 
considering the input load as L h (t)= X x e s '* (which in 

dimensionless form is expressed as i b (t)= (b/mUl 

the response of the system is postulated in the form 

h(l)=H l (s l )X l e v +h.oJ- Substituting h(t) and its 
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derivatives in the governing equation of motion, one 
determines the coefficient of X,e v . 

In a linear aeroelastic system, the system is completely 
characterized by a transfer function H } (si) that contains 
the aerodynamic term as follow: 

)=(**, +rn s , 2 + s ^^bU ^[-isfi/U „] 

d 2 ) 

Herein the Theodorsen's function C, connected with the 
Wagner's indicial function 0(r) via the Laplace's 

transform as c(~is)= s^ (p(r)e~ ST dr ^ has been included 

in the formulation. The terms underscored by the solid 
line correspond to the unsteady aerodynamic loads 
component (circulatory term), while the dotted line 
identifies the terms corresponding to the effect of the 
added mass. When the aerodynamic loads are neglected 
and for s = ico, this result coincides with that of the 
linear FRF, derived via the conventional Modal 
Analysis. 

For purely mechanical systems, in the frequency 
domain, the response analyses via Volterra series have 
been carried out by several authors. In the present study 
an alternative procedure, based on the multivariable 
kernel transforms referred to as Higher-order Transfer 
Functions (HTFs) is pursued. The two above mentioned 
approaches can be correlated each other, and this is 
shown also in this work. Assuming zero initial 
condition, the frequency response functions (FRFs) are 
obtained from the transfer functions TFs, by replacing 
the Laplace transform variable s with jco where co is the 
frequency of the excitation, (Worden et al. 33 ). 

In the present nonlinear aeroelastic system, toward the 
estimation of higher order frequency response functions 
(HO FRFs) that are defined as the multi-dimensional 
Fourier Transform (MDFT) of the Volterra’s kernels, a 
sequence of transfer functions are employed. The 
concept of higher order of FRFs, independent of the 
input to the system, defined from the Volterra series, 
will also be included. 

By the use of the linear frequency-response-function 
the behavior of the linear system is easily 

determined. It will be necessary to find a complete set 
of Volterra kernel transforms H tl (s /f s 2f --s n ) for 

nonlinear systems and for this, in practice, we will use a 
convergent truncated series. 

However, probing the system with a single harmonic 
yields only the information about the value of the 
transfer functions terms on the diagonal line of the plane 
Sj, S 2 , in the Laplace transformed space, where Sj = s 2 > 
However, in order to obtain information elsewhere in 
this space, one should use multi-frequency excitations. 

In the same way, the second order Volterra Kernel can 
be determined applying a load depending on two 


different frequencies expressed as: L h (t)= X,e v + X 2 e SyJ . 

In this case we can express the plunging response in the 
form: 

h(t)= A/, (s, )X ] e' 1 ' + //, (j, )X 2 e‘ 1 ' 

+ H 2 (s i ,s i )X 2 e 2h ' + H 2 (s 2 ,s 2 )* 2 V v 
+ H 2 (j, , j 2 )X , X 2 e (, ' t, ‘' f + H 2 (s, , s, )X 2 X > 

+ h.o.t. .(13) 

Substituting Eq. (13) in Eq. (8) and equating the terms 
containing + ?2 ^ the second order aeroelastic 

Volterra Kernel in the Laplace transformed space is 
obtained: 

H 2 {s\,s 2 ) — —{S\S 2 c } ,2 + &y,2 O'! )^l { S 2 (^1 +£ 2 )’^^) 

where: 

//,(s, +S 2 )=(k hl + (*, + S 2 fm+C lli (s l + J 2 ) 

+ (j, + s 2 )pC La bU_C[- i(s, + j 2 )b/U„ ] 

"2 C s 'i + s 2 ) b ) , (15) 

is the first order Volterra Kernel in the Laplace 
transformed space at the frequency (0|+ Uh. Following 
the same steps, applying the load 
L h (t)= X t e s,r + X 2 e S! ' + X ) e’ , \ equating the terms in the 
form ; remembering that 

//,($, +s 2 )=(^, +( 5 , + s 2 )‘ m + c, (s, + .r 2 ) 

+ (i, + s 2 )pC La bU „C[-i{s, +s 2 )blU„] 

+ ±pb 2 C La (s ] +s 2 ) 2 Y . (16) 

H,(s t +s 2 + .S, )=(<:, + (s, +s 2 +J 3 ) 2 m + C|(j, + s 2 + s 3 ) 

+ ( J i +s 2 +s,)pC La bU„C[-i(s t +s 2 +s,)b/U„] 

+ipt 2 C la (, 1+ 5 2+ 5 3 ) 2 )‘' ,(17) 

the expressions for the third order Volterra Kernel in the 
Laplace transformed space can be cast as: 

2 

H ^{s }t S 2 -,S 7i ) = ~ — {H \ (.?| )// 1 (.V 2 + c Yvb s 2 *3 ) 

+ 2H 2 {s i , s 2 X^Vi 2 C h 2 (* V 1 S 2 V .3 )) 

+ 2(//j ($ 2 )H 2 C? 3 , S t yjc h 2 + C h2 S 2 ( S l + S } )) 

- H , (s 2 )H 2 (s 2 , 5j Xk h2 + c A2 5, (s 2 + j , ))))/ 

("■(*, +J 2 +5,)) .(18) 

Notice that the constants k hl and Cft2 multiply the 
whole expression for and this term vanishes if the 
quadratic term is absent in the aeroelastic governing 
equation of motion. As one of the general properties of 
Volterra’s series, if all nonlinear terms in the equation 
of motion for the system are odd powers of x and y, then 
the associated Volterra series have no even-order 
kernels. As a consequence it will possess no even-order 
TFs. It is also a general property of systems that all 
higher-order TFs can be expressed in terms of H x ■ The 
expressions are function of the system and can be 
obtained using the harmonic probing algorithm. 
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Plunging-Pitching Airfoil 

The governing aeroelastic system of an airfoil 
featuring plunging - twisting coupled motion, exposed 
to a harmonic time dependent external excitation is: 

mh + S a d + ^(c„jh J + k hJ h J )~ L a = L h ’ < 19) 

7=1 

+ +k aj ct i )-M a = 0 - ( 2 °) 

7=1 

Following the steps adopted for 1 DOF, applying a load 
depending on one frequency L h =X l e s ' l > and 
expressing the plunging and pitching displacements in 
terms of transfer functions as: 


a{t) = X , //“ (.V, >'■' + X;H? (s, , .5, > 2 + X, 5 //? (.5, , .5, , .5, y ,(22) 
the relative kernels and the aeroelastic responses can be 
determined. 

The aeroelastic governing system including the blast 
pressure signatures can be expressed in the Laplace 
transformed space as: 

s 2 i+xy«+ 2 CHTr4+(^] I 


2 

+ — 


+ S 2 £ + J"| — 


w 


+— S 2 (f-«d)+— sa=l b {s) ■ 
ft fi 

/ 2 


(23) 


(-0 


. (24) 


(*„ l r a V 2 Z + 2 a + ( 2 C „ +< a/V' 2 

ll 2 / £ a fl ^ l l A t l l 2 . _ 

V ' 1 2 Jrj/i 8r> 

Herein (•)=,£(•), so =,£(£(r)) and d=£(a(/)) 
Following the same steps, applying the loads 
4(0= and L b {t)= X^ 1 ' 1 + X 2 e ht + X 3 e r ’ r > 

equating the terms in the forms X,X 2 ^ {s,+JlJr and 
X l X 2 X i e i ' Sl+Sl+s ' )l the expressions for the second and third 

order Volterra Kernel in the Laplace transformed space 
can be obtained. 


Generalization to Multi Degrees 
of Freedoms Systems (M-DOFs) 

The method shown for 1-DOF and 2-DOF lifting 
surface can be extended to systems featuring multi 
degrees of freedoms, in general, and to a 3-D aircraft 
wing, in particular. 

The method of deriving the n -th order nonlinear 
aeroelastic transfer functions is based upon the fact that 
when the aeroelastic system described by the response 
y(t) (expressed via Volterra series), is excited by a set of 
k unit amplitude exponentials at the arbitrary 
frequencies Sj, S 2 , .... s k , the output will contain 
exponential components of the form: 


>-(')= • (25) 
«=1 m=l 

The presence of nonlinearities causes harmonic 
excitations and sums of harmonics to appear in the 
response of the aeroelastic system. Due to the nonlinear 
formulation, different frequencies can be expected as 
well. 

From the energetic point of view, we can observe that 
// ($ ) produces a single frequency output in response 

to the simple input e *' f . However, because the system is 
nonlinear, H 2 (s lt s 2 ) takes into account the terms that 
produce an output energy corresponding to the sum of 
frequencies cOj+co 2 , or in other words to the input 

Similarly, the third order nonlinear aeroelastic 
kernel, will inject a mix of three input frequencies into 
the total system output. 

RESULTS AND DISCUSSIONS 


To assess the versatility and provide a validation of this 
methodology, a comparison of the predictions of the 
aeroelastic response of nonlinear 2-D lifting surface 
using three approximations are shown in Figs. 2. 



( a ) Time [sec] 



Time [sec] 


Fig. 2 Convergence study involving the first three 
kernels and comparison with the “exact” nonlinear 
aeroelastic response to (a) 1 -COSINE gust pulse and (b) 
to triangular blast load, as shown in the inset 
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The excellent agreement of the predictions, assess both 
the accuracy of the aeroelastic model and also the power 
of the methodology based on the Volterra series and 
indicia! function approach. The first, the second and the 
third approximations of the aeroelastic response to the 
two loads I -COSINE gust load and triangular blast load 
are plotted for different parameters, together with the 
“exact"’ response of the aeroelastic system as obtained 
through digital-computer solution of the nonlinear 
aeroelastic governing equations. Both figures reveal the 
rapid convergence of the approximation. The 
parameters in use for the simulations, unless otherwise 
specified, are chosen as: ( m = \; a = -0.2 ;c h] = 10; 

K =10V OI =10;* OI =10 4 ;c, 2 =10;*,, = 10 7 ;c„ 2 = 10;* n2 =I0 7 ; 
b = \\p=0A25-,U_=0.4U F ;C La =2n)- For the 2-D lifting 
surface encompassing pure plunging, the first three 
aeroelastic kernels in magnitude and phase are depicted 
in Fig. 3 as a function of the frequency, considering that 
(Q = Q)i=(o 2 =(Q i , i.e. the representation is given along 
the diagonal of the plane cu,,^ 2 . As is clearly seen, a 
reduced influence on the response of the third kernel is 
experienced. 


Magnitude [10 A 3] 




Fig. 4 First two aeroelastic kernels of 2-D lifting 
surface, plunging - pitching coupled motions 





Frequency , to ( =<qj=<o>=(o > [Hr) 
Phase fried 



V 



1M lit 


Frequency , to, =0^=035=0^ [Hz] 


Fig. 5 3-D and contour plots of second 
order aeroelastic kernel 



Fig. 3 First three aeroelastic kernels 
of 2-D lifting surface, pure plunging motion 

In Figs. 4 the Volterra’s kernels for the lifting surface 
featuring plunging - pitching coupled motions are 
depicted. Also in this case in the plots include the 
magnitude and phase for the kernels in plunging 

and pitching Hf , in which i identifies the order of 
the kernel. 



Fig. 6 3-D and contour plots of third 
order aeroelastic kernel 
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71 (t) 


7><t) 


A 3-D plots of the magnitude and phase of the second 
order kernel vs. the two frequencies and (0 2 are 
displayed in Fig. 5. The contour plots reveals the 
symmetry of this kernel respect to the diagonal 
represented by =o> 2 . In a 3-D plot, the third order 

Volterra kernel for the case in which &> 3 = CO ] are 
depicted in Figs. 6. 

Response in Time and Frequency Domains 

Determination of subcritical aeroelastic response to any 
time-dependent externally applied load is useful in the 
design of wing structures and of the associated feedback 
control systems. In certain types of nonlinear analysis 
we are only interested in the special case considering of 
T, = t 2 = • • • = r = t 34 - This case can be represented as: 

*(*)■ MV r 2 .-Tj ti=Tj= ... =f<=r - (26) 

This has a corresponding Laplace transform g(s) (so 
called associated transform ) in the single-dimensional 
Laplace transform space: g(s)=£ [g(r)]< The response 
in time can be obtained from H(s ]t s 2 ,-',s fl ) to Find 
G(s) first and evaluate the single dimensional inverse 
Laplace transform g(r). 


This approach is called association of variable 34 . The 
nonlinear aeroelastic response in the time domain is 
depicted in Figs. 7 for a 2-D lifting surface featuring the 
plunging degree of freedom. In this figure the first plot 
represents the linear impulse response that corresponds 
to the convolution integral for the linear analysis. The 
other three plots represent the components of the 
response due to the second and the third order kernels 
and the total response as a combination of the three 
partial responses. The aeroelastic response will be 
presented and validated. The influence of the linear and 
nonlinear stiffness and the damping coefficients on the 
response, not displayed in this paper, reveals that, an 
increase in the damping coefficient contributes to the 
decrease of the response amplitude. 

An increase of the nonlinear damping or of the stiffness 
coefficients contributes to the decrease of the magnitude 
of the kernels and consequently, of the response 
amplitude. This shows that the nonlinearities in the 
stiffness and damping play a beneficial role on the 
subcritical aeroelastic response 15 . 

Figure 8 highlights the effect of the speed parameter V 
(=U/bo> a ) on the lifting surfaces subjected to sonic- 

boom pressure signature as shown in the inset. 


Nonlinear Aeroelastic Responses — Time Domain 


Linear ImpulseResponsewith onekemel Impulse Response Second Order Kernel 




Impulse Respond Third Order Kernel NonttnearlmputoeResponsejCOsy, (0+^ (t)+7j(t) 




Fig, 7 Time-history of the nonlinear aeroelastic response 
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Fig. 8 Influence of the flight speed on the nonlinear 
aeroelastic response to a sonic- boom, as shown in the 
inset (z p = 15sec;r = 2), evaluated with three kernels. 


Herein x p denotes the positive phase duration of the 
pulse measured from the time of impact of the structure; 
r denotes the shock pulse length factor. For r = 2 the re- 
shaped pulse degenerates into a symmetric sonic-boom 
pulse, in the sense that its positive phase has the same 
characteristics as its negative one, and for r = 1 a 
triangular pulse that corresponds to an explosive pulse is 
obtained. It becomes apparent that the amplitude of the 
response time-history (that have been evaluated for 
practical use with three kernels) increases with the 
increase of V . Moreover, in a certain range of speeds, 
as time unfolds, a decay of the amplitude is 
experienced, which reflects the fact that in this case the 
subcritical response is involved. However, for the 
dimensionless speed parameter V greater then the 
flutter speed (this one was determined using the 
linearized aeroelastic system), the response becomes 
unbounded implying that the occurrence of the flutter 
instability is impending. Also in this case the nonlinear 
stiffness and damping coefficients play a beneficial role 
on the subcritical aeroelastic response. 

CONCLUSIONS 

Several issues that concern the nonlinear aeroelastic 
response via Volterra’s series approach have been 
presented. It was also shown that, the method based on 
Volterra series opens large opportunities to approach in 
an unified and efficient way problems of nonlinear 
aeroelastic response and flutter. In addition, following 
the same approach, the character of the instability 
boundary, i.e. benign or catastrophic will also be 
addressed. This analysis will be done by using the 
concept of the first Liapunov quantity as developed by 
Bautin 35 . 


Moreover, this approach can be extended as to include 
also active control capabilities. In spite of this, few of 
these potentialities have been explored yet. 
Comparisons of results carried out via Volterra series in 
conjunction with indicial functions approach and 
classical approach have been provided in Fig. 1 for the 
linearized model. It should also be stressed that 
aerodynamic indicial functions (for incompressible/ 
compressible flow fields 31 ) considered in conjunction 
with Volterra’s series approach can be used as a 
powerful analytical tool for developing unsteady 
aerodynamic models and a unified nonlinear aeroelastic 
model. To the best of the authors' knowledge, with the 
exception of this paper, the problem of the aeroelastic 
response of lifting surfaces to external pulses via 
Volterra’s series and indicial function approach was not 
yet addressed in the literature. 
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